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Abstract. We address a mean- field zero-temperature Ginzburg-Landau, or 4 , 
model subjected to quenched additive noise, which has been used recently as a 
framework for analyzing collective effects induced by diversity. We first make use 
of a self-consistent theory to calculate the phase diagram of the system, predicting 
the onset of an order-disorder critical transition at a critical value a c of the 
quenched noise intensity er, with critical exponents that follow Landau theory of 
thermal phase transitions. We subsequently perform a numerical integration of the 
system's dynamical variables in order to compare the analytical results (valid in 
the thermodynamic limit and associated to the ground state of the global Lyapunov 
potential) with the stationary state of the (finite size) system. In the region of the 
parameter space where metastability is absent (and therefore the stationary state 
coincide with the ground state of the Lyapunov potential), a finite-size scaling analysis 
of the order parameter fluctuations suggests that the magnetic susceptibility diverges 
quadratically in the vicinity of the transition, what constitutes a violation of the 
fluctuation-dissipation relation. We derive an effective Hamiltonian and accordingly 
argue that its functional form does not allow to straightforwardly relate the order 
parameter fluctuations to the linear response of the system, at odds with equilibrium 
theory In the region of the parameter space where the system is susceptible to 
have a large number of metastable states (and therefore the stationary state does 
not necessarily correspond to the ground state of the global Lyapunov potential), we 
numerically find a phase diagram that strongly depends on the initial conditions of the 
dynamical variables. Specifically, for symmetrically distributed initial conditions the 
system evidences a disorder-order transition for a' c < a Cl yielding a reentrant transition 
in the full picture. The location of a' c increases with parameter a and eventually 
coalesces with a c , yielding in this case the disappearance of both transitions. On 
the other hand, for positive-definite initial conditions the order-disorder transition is 
eventually smoothed for large values of a, and no critical behavior is found accordingly. 
At this point we conclude that structural diversity can induce both the creation and 
annihilation of order in a nontrivial way. 
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1. Introduction 

In statistical mechanics, models describing the effect of impurities or heterogeneities in 
the behavior of magnetic systems are gathered under the label of spin glasses [1] when the 
source of heterogeneity affects the local spin interaction (and therefore the interaction 
term in the Hamiltonian takes into account such disorder). Conversely, the so-called 
random field models [1] address those systems where the source of heterogeneity only 
yields an additive heterogeneous term (perturbation) in the Hamiltonian: in this case 
the effect of disorder is akin to subject the system to a random external perturbation. 
In both cases, such sources of heterogeneity typically have slower dynamical evolution 
than the spins (or the dynamical variables), and therefore these sources of randomness 
are said to be quenched. In the last decades a wealth of literature has addressed the 
phenomenology behind spin glasses and random field models, including phase diagrams, 
aging and other dynamical behavior, and comparison with their equilibrium counterparts 
(see [H EJ [3] and references therein). 

In other branches of science the role of disorder in models characterizing the 
dynamical behavior of multicomponent systems has also been addressed in the last 
years. Noticeable examples include the effect that a certain amount of heterogeneity 
in the natural frequencies of Kuramoto oscillators can yield on synchronization [U [5] , 
the paradoxical constructive role that disorder can induce in the formation of ordered 
structures in a plethora of different contexts [HI El El IHl HQl HU H21 HH HH H5l HH1 HZ], 
and the effect that the topology of the underlying network of interactions plays in 
several types of dynamics [TH [191 EDJ [21], to cite some. All these works address similar 
generic questions, namely study the effect of structural disorder in the dynamics of 
multicomponent systems. 

In this paper we will address a paradigmatic example within equilibrium statistical 
mechanics, the Ginzburg-Landau, also called 4 , model [22], in a version subjected 
to such quenched disorder much in the vein of random field models. Although the 
expected role of heterogeneity is that of destroying the ordered state, recent works [7J [9] 
have addressed the positive role of the quenched noise in enhancing the response of this 
model under the presence of an external periodic driving. In [6] the authors studied the 
effects of introducing a quenched multiplicative dichotomous noise, and found that the 
phase diagram is modified and gives rise to the onset of reentrant phase transitions not 
present in the quenched noise free model. 

Here we address the mean-field version of the model subjected to quenched additive 
noise in absence of temperature j9j[23]. First, we present an analytical study of the phase 
diagram by means of a self-consistent theory, both in the non-metastable and metastable 
regions. The theory predicts an order-disorder transition as a function of the quenched 
noise intensity a, with mean field critical exponents equal to those of the thermal 
equilibrium counterpart. We also perform a detailed numerical study of the system 
for different sizes N in terms of finite-size-scaling theory and determine the scaling 
exponents. We show that in the non-metastable region the order parameter fluctuations 
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diverge with an exponent different from the one of the magnetic susceptibility. This 
indicates a violation of the fluctuation-dissipation relation. In order to justify this 
finding, we obtain in closed form an expression for the probability density function of 
the system in terms of an effective Hamiltonian "H e fj(x), and accordingly argue that the 
fluctuations of the order parameter cannot be straightforwardly related to the linear 
response of the system. In the region where metastability does take place, results from 
numerical simulations deviate from the phase diagram found through the self-consistent 
theory and show a strong dependence on the specific initial conditions. Concretely, we 
show that for symmetrical initial conditions, the simulations point out the presence of a 
reentrant phase transition (disorder-order-disorder) with an ordered state whose width 
varies and eventually disappears in the Ising limit, corresponding to a large valued of 
a parameter in the Hamiltonian. This counterintuitive phenomenology supports the 
fact that disorder or heterogeneity can not only induce dynamical disorder but, on the 
contrary, can have an ordering role. Conversely, for positive-definite initial conditions 
the phase transition is smoothed in the same limit, and no critical behavior is found in 
that case. 

The rest of the paper is organized as follows: in section [2] we present the model. In 
section [3] we outline some considerations regarding the presence of metastable states. In 
section H] we derive the mean-field critical exponents associated to the magnetization and 
magnetic susceptibility. In section [5] we numerically study the order-disorder transition 
in the range of parameters where the system lacks metastable states. We provide 
compelling evidence suggesting that the fluctuation-dissipation relation is not satisfied, 
and we argue that a possible reason is that the influence of the average external field 
h on the effective Hamiltonian yielding the probability density function of the system 
cannot be readily stated as % e ff( x ) = %o( x ) + Nmh, being m the magnetization, as it 
happens in equilibrium theory. In section [6] we numerically explore the system's behavior 
in the presence of metastable states and discuss the role of the initial conditions in 
the asymptotic stationary state of the system. We also point out the presence of an 
disorder-order transition induced by diversity in the metastable situation. In section [7] 
we summarize our main results. 

2. Additive Ginzburg-Landau model: preliminary considerations 

We consider a set of N real dynamical variables Xi(t),i = 1, . . . , N whose evolution is 
given by a relaxational gradient flow [21] in a potential V: 

dxi dV{x;rj) 



dt 



dxi 




(1) 
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or, 

Jjl = axi - xf + J2(xj - Xi) + r]i . (2) 

3=1 

The Lyapunov potential V{x;i]) depends, besides on the dynamical variables x = 
(xi, xjv), on a set of variables r\ = (r/ 1; r] N ). Most commonly these variables 
represent white noise of amplitude proportional to the temperature and the model 
defines a class of thermal phase transitions. In this work, however, we take these 
variables to represent quenched noise and the problem then belongs to a class of zero 
temperature random field models. Accordingly, (r]i, . . . ,tjn) are independently drawn 
from a probability distribution g(rj) (which typically will be a Gaussian) of mean h and 
standard deviation a. The model can be thought as describing a set of globally coupled 
heterogeneous units, being a a measure of the amount of diversity or heterogeneity in 
the system. As we are interested in this work in the effect of the diversity, a will be 
taken as a control parameter and we will study the effect that a has on the collective 
properties of the system. 

This model is indeed a discretization of a mean-field version of the well known 
Ginzburg-Landau Hamiltonian for a scalar field x{r) under the presence of a random 
external field 77 (r) [TJ 122] : 

U = J dr(^-^x 2 + ^\Vx\ 2 + ^x 4 -nx^ , (3) 

where, without loss of generality, we have rescaled variables and time such that 
u = 1, C = 1/2. This Hamiltonian provides a coarse-grained description of critical 
phenomena, and its formulation is based on some phenomenological considerations such 
as locality and symmetries (rotational and translational); that is to say, this latter 
expression is not calculated from the microscopic physics, but rather can be understood 
as a coarse-grained description of the magnetization field x. By using the Boltzmann 
weight factor e~ n ^ T , where T is the temperature, this model has been used for instance 
to describe the paramagnetic-ferromagnetic transition (where the Hamiltonian describes 
the coarse-grained magnetization field). In the case of a uniform external field, Landau 
theory elegantly describes a second-order thermal phase transition for this system, with 
mean-field critical exponents ,5 = 1/2, 7 = 1 [H [22] . This Hamiltonian also offers a 
soft-spin description of the Ising model [I]: as a matter of fact, in the limit a — > 00 
one recovers the Ising model (or the Random Field Ising Model (RFIM) in the case of 
having a random external field). In the last decades the RFIM has been extensively 
studied (see PQE5] and references therein), where some specific results include the onset 
of criticality in terms of a second-order phase transition in zero-temperature induced by 
the disorder of the random field, with mean-field critical exponents [261 12Z] as in the 
thermal counterpart [28]. Several other features such as hysteresis, avalanche dynamics, 
or return point memory effects, to cite a few, have been studied within the RFIM, both in 
analytical (renormalization-group) and numerical (finite-size scaling) terms [251129]. The 
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properties of the model have also been studied in the context of domain growth dynamics 
both in the Ising limit |30[ |3T| [32] or using the full Ginzburg-Landau Hamiltonian (33] . 

3. On the presence of metastability 

From the dynamical point of view, it has already been said that the evolution is 
relaxational in the Lyapunov potential V. Hence, the absolute minimum (or ground 
state) of V located at x e . . . , x^) must be considered as the global attractor of 
the dynamics. It is obvious that the value of x will depend on the specific realization 
of the quenched noise variables (771, . . . ,n N ). On the other hand, the solutions of the 
differential equations (j2|) tend to values xf = \imt-too %i{t) which might or might not 
coincide with X{. If the potential V has a single minimum, then the dynamics always 
leads to x, but if there are additional, metastable, minima, then the asymptotic solution 
x st depends on the initial condition x(t = 0) as it might get stuck in one of them. The 
presence and relevance of these metastable minima depends in general (and besides the 
particular realization of the quenched-noise variables) on the value of the parameter a 
and the number of variables N. 

In order to find the absolute minimum x one needs to solve the system of N coupled 
algebraic equations: 

1 N 

Q = axi-Xi+ — y]{xj - Xj) + ijj. (4) 

3=1 



The solution is greatly simplified if one introduces the magnetization m as 

N 

1 ■ 

m 



]yX>> ( 5 ) 

i=i 

and then writes Eq.(j4|) as: 



m + r]i = (1 — a)xi + x\ . (6) 

This equation allows one to find X{ as a function of m and T]i (in fact as a function of 
m + rji). The explicit solution, X{ = x(m + 77^) can be replaced in the definition of the 
magnetization to obtain a self-consistency equation: 

1 - 

m = —^x(m + r]i). (7) 

i=l 

The problem has been reduced from the simultaneous solution of the iV coupled 
equations (@J, to the solution of a single one © although, in general, all possible 
solutions 777^,777^,... of this equation have to be obtained numerically. For a 
given solution one can then find the respective values of xf 1 ^ using the function 
x\ n ^ = x(m( n ) +T)i). In order to analyze the structure of the possible solutions of Eq.flT]), 
it is convenient to split the discussion in the cases a < 1 and a > 1. 
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3.1. Case a < 1 

This is the simplest case. A graphical analysis shows that Eq.([6]) has a unique real 
solution x-i = x{m + rji) (see Appendix). Even in this case, it is possible that Eq.([7J has 
more that one solution for m. This is typically the case for small values of N. See an 
example in Fig. [U 




Figure 1. Lyapunov potential V(xi,X2) as defined in Eq. (TT]) for N = 2, rjx = 
—0.48, r\2 — 0.5 in the cases a = — 1 (left), a = 0.8 (center) and a = 2.8 (right). While 
the case a = — 1 displays a single minimum, in the case a = 0.8 there are 3 minima 
(2 metastable) and 2 maxima, whereas for a — 2.8 there are 4 minima (3 metastable) 
and 5 maxima. 



However, as N increases the number of metastable solutions decreases greatly. In 
fact, it is possible to prove that in the thermodynamic limit, N — > oo, Eq.([7J can have 
only either one or three solutions depending on the values of a, h, a. The proof replaces 
the sum over N by an integral over the probability distribution of the quenched-noise 
variables: 

m = J drj g(j])x(m + rf). (8) 

Let us assume that the probability distribution g(j]) has a generic form g(rj) = • 
Henceforth, all numerical results will use the Gaussian distribution G(z) = ~^ e ~ z2 ^ 2 - 
A change of variables leads to: 

77i= / ' d£G(Z)x(m + h + a£) = F a {m + h). (9) 

As F a (z) is a monotonously increasing function satisfying F a (0) = and with a 
sigmoidal shape [36], there will be only one solution for m for all values of h if 
the derivative satisfies FL(0) < 1. On the other hand, for FL(0) > 1 there will 
be either one or three solutions depending on the value of h. This analysis mimics 
that of the Weiss mean-field theory [37] and allows one to compute the magnetization 
m(h; a, a) as a function of the mean value of the disorder h and the parameters a, a. 
It displays usual critical phenomena and hysteresis. The critical point is defined by 
the condition F' a (§) = 1 and can be achieved by varying a or a. It is possible to 
show that F^ =0 (0) = 1/(1 — a) and, since F^(0) is a decreasing function of a, the 
condition F^(0) = 1 can never be achieved for a < 0. This was a priori obvious since 
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in that case the Lyapunov potential in the absence of quenched noise has the global 
minimum at x^ = 0, Vi, already a disordered state. Some numerical values (for the 
Gaussian distribution) for the location of the critical diversity function of a 

are: (a = 0.1, 07. = 0.19616), (a = 0.5, a c = 0.50041), (a = 2/3, a c = 0.595233). In 
the case a = 1, the Cardano formula simplifies to x = (m + h) 1 ^ 3 and it is possible to 
perform analytically the integrals (again for a Gaussian distribution for the quenched- 

noise variables) with the result [9] (a = 1, o~ c 



rq/6) 

2 1 /3 37r l/2 



3/2 



0.7573428 . . .) 



3.2. Case a > 1 

The problem in this case is that the cubic equation ([6]) can have either one or three 
real solutions depending on whether the discriminant Aj = 27(m + rji) 2 + 4(1 — a) 3 is, 
respectively, positive or negative. Besides, as before, several values of m can satisfy the 
self-consistency Eq.([7]). When there are three solutions for Xi, (Aj < 0, this requires 
a > 1) it is not clear a priori which one to chose in order to substitute in the self- 
consistency relation ([71). A possibility is to compute the Lyapunov potential V for each 
of the possible solutions. However, since the maximum number of solutions can be 
as large as 3^, this is not possible to carry out in practice for large N. The answer 
arises when one realizes that the dynamical equation for Xi(t) can be written also as 
relaxational in a local potential v i(xt, m)\ 

dxi dvi(xi,m) 
dt dxi 

Vi = — ^— x t + -Xi - [m + r)i)Xi + — . (10) 

The solutions x(m + are nothing but the extrema of this local potential. Now we 
notice that the Lyapunov potential can be written as sum of the local potentials: 

N 

V{xx, . . .,x N ) = yvj{x h m). (11) 

1=1 

Therefore, the absolute minimum of V is achieved by choosing in each case the solution 
x(m + rji) that minimizes the local potential Vi(xi,m). Explicit expressions for the 
function x are obtained using Cardano's formula and are given in the Appendix. 

The process to find the absolute minimum x of the Lyapunov potential proceeds, 
as before, by finding first m after solving numerically the self-consistency equation ([7J, 
but using the correct function x(m + i]). Similarly, the integral equation ([8]) can be used 
to find the magnetization m(h; a, a) in the thermodynamic limit. The phenomenology 
of the solutions is similar to what was found in the case a < 1 and will not be repeated 
here. 

An important difference, however, with the case a < 1 is that now the Lyapunov 
potential displays a large number of metastable minima for all values of iV and, 
consequently, also in the thermodynamic limit (a recent study for the metastable states 
of the zero-temperature RFIM has been carried on in [M, 35J). Therefore, starting from 
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arbitrary initial conditions, the asymptotic solution of the evolution equations xf will 
in general differ from the values Xj of the absolute minimum. It will be shown that new 
phase transitions occur when one looks at the magnetization values that derive from 
the stationary solution. 



4. Critical behavior 



We have seen that this mean-field model displays a second order phase transition 
between an ordered state {\m\ > 0) and a disordered state (m = 0) at a critical 
value of the diversity a c . In order to derive the critical exponents of such transition, 
we consider the self-consistency Eq. ((9j) and expand F a (m + h) in a Taylor series. 
Since x(—m — h) = —xim + h) (see Appendix) and assuming that the distribution 
of noises is symmetric with respect to the mean value, G(— £) = G(£), the function F a 
is antisymmetric F a (—m — h) = —F a (in + h) and we get: 

m = ai(a)(m + h) + a 3 (a)(m + h) 3 + . . . (12) 

with afc(cr) = (0)/k\. Hence, the magnetization at h — is: 

for a > a c , 



\m\ 



ai (a) (13) 
tor a < a r . 



03(0) 

As i^(0) — 1 changes sign at a = a c , we can expand ai(a) = 1 + ai(a c — a) + — 
Accordingly, close to the transition the spontaneous magnetization behaves as \m\ ~ 
(<t c — a)?, with a critical exponent (3 = 1/2, as in Landau's treatment of the thermal 
phase transition. 

Dm 

To compute the critical behavior of the susceptibility Xh = ttt , we take 

dh h=0 

the derivative of both sides of Eq. (1121) and set h = 0. This leads to Xh = 
ai(a) + 3a 3 (a)m 2 

— — — — — -. Replacing Eq.( T131) and aAa) = 1 + aAa c — a) + . . . we find 

1 — ai(cr) — 6a^{a)m z 

the critical behavior: 

Xh = A±\a -a^ 1 (14) 

with critical amplitudes = l/(2«i) for a < o c and A + = l/ot\ for a > o c . Therefore 
the susceptibility critical exponent is 7 = 1, the same, not surprisingly, than in Landau's 
theory. 



5. Numerical results for a < 1: violation of the fluctuation-dissipation 
relation 

In this section we present the results coming from the numerical integrations [38] of 
the dynamical equations (J2j) in the case a < 1. The objective is twofold. First, by 
comparing with the analytical results valid in the thermodynamic limit, we want to 
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a=2/3 a=2/3 
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0.25 0.5 0.75 I 0.56 0.6 0.64 0.68 0.72 



Figure 2. Left panel: Order parameter mo as a function of the diversity a for a = 2/3. 
The symbols correspond to the numerical integration of the dynamical equations ([2]) 
for different system sizes N and a Gaussian distribution (zero mean, standard deviation 
a) of the quenched noises. The solid line is the magnetization m obtained by solving 
the self-consistency Eq. {[9j) for h = 0. Right panel: Order parameter fluctuations, x, as 
a function of the diversity a, for the same system sizes as the left panel (the vertical 
axis is in logscale for presentation purposes). 

check the importance of the metastable states that appear for finite N. Second, we will 
use the theory of finite-size scaling in order to determine the exponents of the transition. 
We will show that there is a violation of the fluctuation-dissipation relation in the sense 
that the magnetic susceptibility can not be computed as the ensemble fluctuations of 
the magnetization. By ensemble average ((■ ■ •)) we mean an average with respect to 
realizations of the random quenched-noise variables as well as with respect to the initial 
conditions x(£ = 0). However, for the range of values of system size N employed in the 
simulations, N > 10 3 , there was hardly any dependence on the initial condition for a 
given realization of the random variables. This shows that metastable states either do 
not exist or it is rare to get trapped in them for this range of values of a and N. In the 
left panel of Fig. [2] we plot the order parameter mo as a function of the diversity a for 
the value a = 2/3. As usual [39], the order parameter is defined as the ensemble average 
of the absolute value of the magnetization m = {{\m\)) computed from the stationary 
values as m = ■h^i^xf. As predicted by the self-consistent treatment explained 
in previous sections, there is a phase transition from an ordered (ferromagnetic-like, 
m > 0) to a disordered (paramagnetic-like, m = 0) phase as a function of o. The 
transition is smeared out by finite-size effects, but it approaches the solution of the 
thermodynamic limit and the transition point a c as the system size N increases. In the 
right panel of this figure we plot the normalized fluctuations of the order parameter, 
defined as x = l(( m2 )) ~ (d m l)) 2 ] as a function of the diversity a. These fluctuations 
have a maximum in the neighborhood of a c and, as shown in the right panel of Fig. |3l 
they increase with increasing N as xi a c) ~ N b with b ~ 2/3 for different values of the 
parameter a, and hence diverge in the thermodynamic limit. As shown in the left panel 
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Figure 3. Log- log plots of the order parameter mo (left panel) and the fluctuations 
X (right panel) as a function of system size N for different values of a, at the 
corresponding critical point a c (a). In all cases we find a good fit to a power-law 
behavior: m ~ N~ c and x ~ N b with c = 0.16 ± 0.01 and b = 0.66 ± 0.02. 



a=2/3 
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n=3-io-; 

N=10* 
K=3-10, 
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N=10^ 
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-10 -5 5 
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Figure 4. Data collapse of the order parameter mo (left panel) and the fluctuations x 
(right panel) according to the finite-size scaling relations mo(cr, N) = N~ v / 2 f m {N v (1 — 
a/a c )) and x(<r,N) = N b f x {N v {l - a/a c )) using v = 1/3, b = 2/3. The goodness of 
the collapse is an evidence supporting the validity of the scaling relations. 



of the same figure, the order parameter at the critical point decreases as m (cr c ) ~ N~ c 
with c m 1/6 and tends to zero in the thermodynamic limit. 

Data for a range of values around the critical region can be collapsed through 
standard finite-size analysis (101 E] according to the scaling laws: m (a,N) = 
N- c f m (N v (l - a/o-c)) and x(cx,iV) = N b f x (N v (l - cr/a c )) with appropriate scaling 
functions f m and f x . A good fit, see Fig. HJ is obtained with v = 2c ~ 1/3. Note 
that this scaling form implies that in the infinite-size limit mo (a) ~ \o~ — a c \^ and 
x(c) ~ \o~ — o- c |~ 7 , with critical exponents = c/v = 1/2 and 7 = b/v = 2. We have 
also performed a finite-size scaling of the fluctuations of the stationary value of the 
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energy (global potential) at the critical disorder a c (a = 2/3) = 0.595233, according to 
which one finds a value for the critical exponent of those fluctuations a ~ 0, the same 
than the (thermal) mean-field result for the specific heat (data not shown). 

While the result of the previous section proved that the susceptibility Xh has a 
critical exponent 7 = 1, the numerical simulations suggest that the fluctuations \ 
diverge close to the critical point as a power law with a different exponent 7 = 2. 
This seems to constitute a violation of the fluctuation-dissipation relation. Since we 
have restricted this analysis to the range a < 1, this violation does not seem to be 
related to typical situations of metastability, absence of time translation symmetry or 
aging [H El [3]. Furthermore, the hyperscaling relation 2/3 + 7 = d c u, that holds in the 
mean-field regime or for d > d c , is satisfied using 7 = 2 as it is known |32] that the 
upper critical dimension is d c = 6 and v = 1/2. 

To explain this discrepancy, we note that the fluctuation-dissipation relation is 
obtained typically for a system in the canonical ensemble at temperature T and whose 
probability density function (pdf) is / x = Z~ x exp(— H/T), with a partition function 
Z = J dx exp(— %(x)/T), being H, the Hamiltonian of the system. If the Hamiltonian 
contains a magnetic interaction "H(x) = "Ho( x ) + Nmh, one can prove the fluctuation- 
dissipation relation between the magnetic susceptibility Xh and the fluctuations of the 
magnetization (m): 



where (■ ■ •) denotes an average with respect to the probability distribution / x (x). 

In our case, there are two averages: with respect to initial conditions and with 
respect to realizations of the random variables rj. We have already argued that for 
a < 1 and large values of N, the results are largely independent of initial conditions, so 
all that contributes to the ensemble average ((••■)) are the noise variables. As there is a 
one to one correspondence between the stationary values x and 77 we can write the pdf 
of x in terms of the pdf of 77: 



If we take the r\i 's to be independently distributed Gaussian variables, we have 




(15) 




(16) 



N r 



1 



> vn) = n 



exp(-(7 7i - h) 2 /2a 2 ) 



(17) 



i=i L 



As Eq.flSJ) implies 





3=1 



it is possible to compute the determinant of the Jacobian matrix J;.,- = 



\ j=l J / i=l 



(19) 
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Figure 5. Numerical results of the average magnetization as a function of diversity 
a, for a system of N — 16384 coupled variables for different values of a > 1 (for the 
numerical integration of Eq(2j initial conditions are drawn from a symmetrical uniform 
distribution U[—S, +5]). Note that depending on the specific value of parameter a, three 
different behaviors take place: (/) an order-disorder transition at a c for a = 1, 1.4, (J/) 
a reentrant phase transition formed by a disorder-order transition at a' c coupled to an 
order-disorder one at a c for intermediate values of a = 1.8, and the absence of 

any transition to an ordered state for the larger value a = 2.4. 



Replacing Eqs.f U7tfT9|) in Eq.f lT6|) . one can write the pdf of x as the exponential of an 
effective Hamiltonian /x(x) = Z~ x exp(— H e s), with: 
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3 =i * 



+ 



N 
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i=i L 



[(1 — a)xi + — m — h\" 
2^2 



- In (3x 2 i + 1 - a) 



(20) 



However, as it can not be splitted in the form H e e. — Ho + Nhm, it is not possible (at 
least in a trivial manner) to relate the susceptibility to the fluctuations of the order 
parameter. 



6. Numerical results for a > 1: dependence on the initial conditions 

In the case a > 1 the presence of metastable states is relevant as the dynamics usually 
gets trapped in one of them. Therefore, in general, the asymptotic values x st depend on 
initial conditions and the absolute minimum of the potential V might not be reached. 
Accordingly, deviations from the self-consistent theory are expected to appear. In 
this section we will study this case and show that a new phenomenology can appear 
depending on the particular value of a and the distribution of the initial condition 
x(i = 0). For the sake of concreteness, we have focused on two types of initial conditions: 
symmetrical and positive-definite. 
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Figure 6. Left panel: Phase diagram of the system, where the symbols correspond 
to the values of critical points a c (associated to the order-disorder transition) and 
a' c (associated to the disorder-order transition) as a function of a, for a system of 
N = 16384 (derived numerically integrating Eq(2] with initial conditions drawn for a 
symmetrical uniform distribution [—6, +5]). In the region a > the system evidences 
an order-disorder phase transition at <r c , the location of this transition increasing with 
a. The values of a c (in the thermodynamic limit) can be derived from the self-consistent 
theory as those satisfying F' a (0) = 1, and are represented by the solid line. In the 
region a > 1 the system presents metastable states even in the thermodynamic limit 
and the solid line refers to the location of phase transition derived from the analysis of 
the ground state of the Lyapunov potential. At odds with the self-consistent theory, we 
numerically find for intermediate values of a the coexistence of two phase transitions 
(reentrant transition) , where the location of both critical points converge for increasing 
values of a until coalescence. At this point the ordered state is completely destroyed for 
all values of a. Right panel: Same diagram as for the right panel, when the numerical 
integration of Eq|5]is performed with initial conditions drawn for a uniform distribution 
in the positive-definite interval [0, 25]. In this case, the phase transitions disappear for 
a > 1.4 as in this case the order parameter mo tends to zero smoothly with a, see right 
panel of FigjSl 

6.1. Symmetrical initial conditions 

The initial values Xi(t = 0),i = 1,...,N, are independently drawn from a uniform 
distribution in the interval [—5, +5], for a given value of 5. In Fig. El we plot the average 
magnetization m = ((|^|)) as a function of diversity a for different values of a and 
system size N = 16384 for 5 = 2.5. The data have been averaged over 10 2 initial 
conditions for x(£ = 0) and then over 10 2 realizations of the quenched noise variables 
(10 4 averages in total). At variance with the case a < 1 (which is also shown in the 
figure for comparison) we find three possible scenarios: (i) for a > 1 (weak metastable 
regime, a = 1.4 in the figure) one observes the same phenomenology as for a < 1: 
an order- disorder transition at a critical value a c (a). (ii) For larger values of a, the 
former transition is still present at a c , but a new transition (from a disordered state 
m = to an ordered one m Q > as o increases) is found at o' c < a c , see the curve 
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Figure 7. Le/t panel: Numerical results of the order parameter as a function of a, for 
different system's size and a — 1.8, where a reentrant phase transition takes place (for 
the numerical integration of EqlJJ initial conditions are drawn from a symmetrical 
uniform distribution U[—8,+S]). Exact results from the self-consistent theory are 
represented in the solid line. The deviations from the theory are related to the fact 
that the system does not reach the ground state of the Lyapunov potential as it gets 
trapped in metastable states. Right panel: Fluctuations of the order parameter as 
a function of a for the same system as the right panel. Fluctuations have a peaked 
maximum that scales with system's size close to both transition points. 



corresponding to a = 1.8 in the figure. In this case, we find the counterintuitive result 
that a certain level of diversity in the quenched noise is needed to induce order at 
a = a' c , whereas a large level of diversity destroys again the ordered state (reentrant 
phase transition), (iii) Finally, for increasing a, a' c increases and a c decreases, eventually 
coalescing for a > a c ~ 2.4, where the ordered state disappears. Thus, for large values 
of a, the system does not evidence any transition and the stationary phase is always 
the disordered one. We point out that in the curve for a = 2.4, the magnetization is 
not exactly zero for intermediate values of the diversity due to a finite-size effect: m 
decreases and approaches zero for all values of a as the system size increases, something 
that does not occur in cases (i) and (ii). All these features are illustrated in the phase 
diagram plotted in the left panel of Fig. [6j (i) For 1 < a < 1.6 the usual order- disorder 
transition appears, although the value of a c is smaller that the one derived from the 
analysis based upon the structure of the global attractor x. (ii) For 1.6 < a < 2.4 there 
is a new transition from a disordered to an ordered state at a value a' c < o c . (iii) Finally, 
for a > 2.4 the only phase encountered is the disordered one. 

In order to characterize the transitions that occur in region (ii), we have run 
extensive simulations for different system sizes in the case a = 1.8. The order parameter 
m is displayed in the left pane of Fig. [3 By looking at the difference with the 
magnetization curve derived from the theoretical analysis, it is clear from this figure that 
the system is not able to reach the absolute minimum neither for small or large diversity 
a. We observe at both transitions the same qualitative dependence with system size 
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that was discussed in the case a < 1. As we don't have now a theoretical prediction for 
o' c or o c the numerical analysis of the data is much less conclusive. Pseudo-critical points 
o~ c (N) and o~' c (N) can be defined as the location of the maximum of the fluctuations x 
of the order parameter, see the right panel of Fig. [7J The fluctuations scale roughly as 
x(cr' c (N)) ~ N b ' and xMN)) ~ N b with b' « b « 0.9. However, it is difficult to obtain 
reasonably good quality fits of the data to the standard finite-size-scaling relations used 
in the case a < 1. Furthermore, the data show a dependence on 5 (data not shown) 
such that a c and cr£, adopt different values for small 5 but saturate for 5 > 2.5. 

Summing up: if the initial conditions are distributed in a symmetrical interval, the 
order region is much reduced with respect to the predictions based upon the structure 
of the ground state. There is a region in parameter space where the system undergoes 
what appear to be well defined phase transitions, from disorder to order and back to 
disorder at a' c and a c , respectively. The order- disorder transition (a c ) is related to the 
one found in the regime a < 1, while the disorder-order transition (at o' c < a c ) is a new 
behavior whose nature is genuinely metastable. For a > 2.4 the system is never in the 
ordered state. 

6.2. Positive- definite initial conditions 

The initial values Xi(t = 0),i = 1,...,N, are independently drawn from a uniform 
distribution in the interval [0,25], for a given value of 5. Obviously, by symmetry 
reasons, the same results would be obtained in the initial conditions were drawn from 
the interval [—25,0]. In Fig. [8] we plot the average magnetization mo = ((|^|)) as a 
function of diversity a for different values of a = 1.2 (left panel) and a = 1.8 (right 
panel), for different system sizes N and values of 5. These two values of a evidence 
slight different behaviors: for a = 1.2, while the sharpening finite-size effect of the 
magnetization is hardly seen in the plot, the fluctuations still increase with system size 
close to the transition (data not shown), what suggests the presence of a phase transition 
in the thermodynamic limit. Note that the dependence on the width of the initial 
condition 5 is very weak and results are basically indistinguishable for 5 > 0.5. On the 
other hand, for a = 1.8 there is hardly any dependence on the system size both for the 
magnetization and its fluctuations. The magnetization mo tends to zero smoothly with 
a and the fluctuations do not increase with system size (data not shown): the transition 
is smoothed and no critical behavior is present. Again, there is a dependence with the 
value of 5 for small 5 but the curves for 8 = 2.5 and 8 = 5.0 are indistinguishable from 
each other. Summing up, for positive-definite initial conditions, the phase transition 
from order to disorder disappears at a value a ~ 1.6 (the actual value depends of the 
width 5), such that the system shows always some degree of order for a > 1.6 (see the 
right panel of Figj6]). In this sense, the ordered region is enhanced with respect to the 
predictions based upon the structure of the ground state. 
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Figure 8. Numerical results of the order parameter as a function of a, for different 
system's size N — 4096,8192,16384, and a value a — 1.2 (left panel) and a = 1.8 
(right panel). For the numerical integration of Eq|2j initial conditions are drawn from 
a positive-definite uniform distribution f7[0,2<5], with S — 0.1,0.5,2.5,3.0. The effect 
of the interval size saturates for approximately 6 > 0.5 and 2.5 for the left and right 
panel respectively. While finite size effects in the magnetization are hardly observed 
for a = 1.2, fluctuations still increase with system size close to the transition. On the 
other hand, for a — 1.8 no finite-size effects are observed, neither for the magnetization 
nor for its fluctuations: the transition is smoothed and no critical behavior is observed. 



7. Conclusions 

In this work we have studied the mean-field version of a Ginzburg-Landau, or 4 , model 
with additive quenched noise at zero-temperature. The model, that has recently been 
proposed in the framework of collective behavior induced by diversity [7J |9] , is a field 
version of the random field Ising model well studied in the literature. As a function of 
diversity a, a self-consistent theory predicts the presence of an order-disorder transition 
at a critical value a c , with mean field critical exponents that are equal than those of 
Landau's theory of thermal phase transitions. Numerical integrations of the dynamical 
equations (j2J) are also performed for comparison. In the range of parameters where the 
system lacks metastable states (a < 1), finite-size scaling relations show that the order 
parameter fluctuations diverge quadratically, rather than with 7 = 1 as in thermal, 
equilibrium, phase transitions. This suggests a violation of the fluctuation dissipation 
not associated to metastable effects such as lack of time translational invariance or aging 
[H El [3]. To explain this fact, we compute an effective Hamiltonian and argue that it 
cannot be readily expressed as "H e g = "Ho + Nhm: as a consequence, the fluctuations 
of the order parameter cannot be straightforwardly related to the linear response, as 
it happens in equilibrium theory. In the range of parameters where metastability is 
likely to appear (a > 1), stationary values typically do not reach the minimum of the 
Lyapunov potential, and accordingly numerical results deviate from the self-consistent 
theory, showing a strong dependence in the initial conditions. For a symmetrical 
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distributed initial condition in the interval [—5, +5] , the ordered region is much reduced 
with respect to the predictions based upon the structure of the ground state of the 
potential. Furthermore, there is a region of values of a for which a new transition from a 
disordered to an ordered state takes place at a' c < a c . In this case, diversity can not only 
destroy an ordered state but also induce order from a disordered metastable state. This 
new transition is genuinely metastable, and its location increases for increasing values 
of a, until coalescing with a c , where the ordered phase completely disappears. On the 
other hand, when the initial condition is distributed in [0, 25], large enough values of a 
destroy the critical behavior of the order-disorder transition and some degree of order 
remains at every value of the diversity a. 

We conclude that structural diversity can induce both the creation and annihilation 
of order in a nontrivial way, and deeply modify the dynamics of the diversity-free system 
counterpart. On the other hand, the apparent violation of the fluctuation-dissipation 
relation should be further investigated; at this point we can conclude that to directly 
relate the order parameter fluctuations to the linear response of a system can be tricky, 
even in the absence of metastability. This is particularly relevant in problems involving 
the estimation of critical exponents in nonequilibrium phase transitions. 
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Appendix: solutions of the cubic equation 

We give explicit expressions for the function x{m + r/) defined as the convenient real 
solution of the cubic equation ax + x 3 = z, where a = 1 — a and z = m + r\. 

In the case a > there is only one real solution to this equation as given by 
Cardano's formula 



x(z) = u — a/(3u), 




For a < 0, the same formula applies if the discriminant A = 27 z 2 + 4a 3 is positive 
A > 0, i.e. z (£ (-2(-a/3) 3/2 , +2(-a/3) 3/2 ) . Otherwise, out of the three real solutions, 
the one that minimizes the local potential v(x) = ^x 2 + \x A — zx is obtained using the 
trigonometric form of Cardano's formula: 

x(z) = 2 sign(2) \J cos ^-arccos 

where the arccos function takes values in the principal branch [0,7r/2] of. Note that, in 
every case, the function x is antisymmetric x(z) = —x(—z). 
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